o 
o 

(N 



> 
oo 



On second order shape optimization methods for electrical impedance 

tomography. 



L. Afraites, M. Dambrine and D. Kateb* 



February 1, 2008 

in 



Abstract 



u 
o 

This paper is devoted to the analysis of a second order method for recovering the a priori un- 
known shape of an inclusion uj inside a body Q from boundary measurement. This inverse problem 
- known as electrical impedance tomography - has many important practical applications and hence 
has focussed much attention during the last years. However, to our best knowledge, no work has yet 
considered a second order approach for this problem. This paper aims to fill that void: we investigate 
the existence of second order derivative of the state u with respect to perturbations of the shape of 
the interface dui, then we choose a cost function in order to recover the geometry of do; and derive 
the expression of the derivatives needed to implement the corresponding Newton method. We then 
investigate the stability of the process and explain why this inverse problem is severely ill-posed by 
proving the compactness of the Hessian at the global minimizer. 
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1 Introduction and statement of the results. 

- i—i ■ 

Let be a bounded open set with smooth boundary in M 2 or M 3 . Consider a L°° function a such 



that there exists a real c with u(x) > c > 0. Consider the elliptic equation 

- div (a(x)Vu) = in CI, (1) 
with the Dirichlet boundary condition 

u = f on Oil. (2) 

Define the Dirichlet-to-Neumann map as 



Act : / a (d a u] 



| an 



where u solves ((T]),([2|) and n is the outer unit normal vector to dfl. The inverse conductivity problem 
of Calderon is to determine a from A ff . Electrical impedance tomography aims to form an image 
of the conductivity distribution a from the knowledge of A CT . When a is smooth enough, one can 
reconstruct a from A CT (see the works of Sylvester and Uhlmann [3T] , Nachmann [T51 US] and Novikov 
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P~7]). When the conductivity distribution is only L°°, Astala and Paivarinta have recently shown in 
[3] that, in dimension two, the map A ff determines a G L°°(f2). 

We are interested in a particular case of that problem: when a body is inserted inside a given 
object with a distinct conductivity, the question of determining its shape from boundary measurement 
arises in many fields of modern technology. In the context of the inverse problem of conductivity 
of Calderon, we restrict the range of admissible conductivity distributions to the family of piecewise 
constant functions which take only two distinct values <ti,02 > which are assumed to be known. 
The conductivity distribution is then defined by an open subset lu as 

o" = oiXn\w + cr 2 Xw (3) 

Here, the only unknown of the problem is lu a subdomain of f2 with a smooth boundary duu; its outer 
unit normal vector is denoted by n. The notation Xuj (respectively. Xfi\w) denotes the characteristic 
function of lu (respectively, fl \ lu). The second main difference arises from practical considerations: 
it is unrealistic from the point of view of applications to know the full graph of Dirichlet-to-Neumann. 
Therefore, we will assume that one has access to a single point in that graph. This non destructive 
testing problem is usually written from a numerical point of view as the minimization of a cost 
function: typically a least-square matching criterion. Many authors have investigated the steepest 
descent method for this problem [T31 [7J [TOl HH H] with the methods of shape optimization since the 
unknown parameter is a geometrical domain. 

This work is devoted to the study of second order methods for this problem that has only be 
considered before for simplified models in [SE]- By introducing second order methods, one aims to 
reach two distinct objectives. 

• On one hand, we provide all the needed material to design a Newton algorithm. We will give 
differentiability results for the state function and for the objective that we have chosen to study 
in this work. Nevertheless, we point out that the discretization of a Newton method for this 
problem turns out to be very delicate; this is why, in the present paper, we will neither discuss 
about this problem nor present numerical examples. This topic is actually the main objective 
of a work in progress. 

• On the other hand, we analyze rigorously the well-posedness of the optimization method. This is 
justified by the huge numerical literature devoted to the numerical study of this question in the 
field of inverse problems; the numerical experiments insist on the ill-posedness of this problem. 
We will explain the instability in the continuous settings in terms of shape optimization. We 
show that the shape Hessian is not coercive -in fact its Riesz operator is compact - and this 
explains the unstability of the minimization process. 

Let us describe the precise problem under consideration and the notations. We consider a bounded 
domain Q, C M. d (d = 2 or 3) with a C 2 boundary. It is filled with a material whose conductivity is a\ 
and with an unknown inclusion lu in £1 of conductivity oi ^ a±. We search to reconstruct the shape 
of lu by measuring on dCl, the input voltage and the corresponding output current. In the sequel, we 
fix do > and consider inclusions lu such that lu CC O,d = {x (z lu, d(x, dQ) > do}. We also assume 
that the boundary dui is of class C A a . The inverse problem arises when one has access to the normal 
vector derivative of the potential u that solves (U)-© when the conductivity distribution is defined 
by ([3|) . Assume that ones knows 

aid n u = g on (4) 

then the problem (JTJ)- (J2J>- (01) is overdetermined. The electrical impedance tomography problem we 
consider is to recover the shape of lu from the knowledge of the single Cauchy pair (f,g). 

In order to recover the shape of the inclusion lu, an usual strategy is to minimize a cost function. 
Many choices are possible; however it turns out that a Kohn and Vogelius type objective leads to a 
minimization problem with nicer properties than the least squares fitting approaches (we refer to [1] 
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for a comparison of different objectives with order one methods and to [5] for the case of a perfectly 
insulated inclusion). Therefore, we study such a cost function in this work. 

Let us define this criterion. Its distinctive feature is to involve two state functions Ud and u n : the 
state Ud solves CU)- d2J) while u n solves (JTJ)-(|4j). The Kohn -Vogelius objective Jkv is then defined as: 



Jkv{u) 



(T\V{u d - u n )f 



(5) 



Let us sum up the results of this paper concerning the minimization of this objective. We first prove 
differentiability results for the state Ud- In the sequel, we use the convention that a bold character 
denotes a vector. If h denotes a deformation field, it can be written as h = h T + h n n on du>. Note 
also that in the following lines, n denotes the outer normal field to doj pointing into Cl\ZJ. Hence, for 
x G duo, we define, when the limit exists, w ± (a;) (resp. (d n u) ± (x)) as the limit of u(x ± tn(x)) (resp. 
(Vu(x ± tn(x) , n(x))) when t > tends to 0. Note that h r is a vector while h n is a scalar quantity. 

The admissible deformation fields have to preserve dd and the regularity of the boundaries: 
therefore the space of admissible fields is 



W = {heC 4 ' Q ( 



J ),Szw(h)cf! do }. 



The following result concerns the first order derivative of the state functions Ud and u n . It was 
derived in [3 EH [Q . 

Theorem 1 Let Vl be an open smooth subset ofM. d (d = 2 or 3) and let lo be an element of ' dd with 
a boundary of class C 4,a . Then the state functions Ud and u n are shape differentiate; furthermore 
their shape derivative u' d and u' n belongs to H 1 (S^ \ ZI7) utffu) and satisfy 



[ad, 



in SI \ to and in lo, 
h n — d n u~, on doj, 

[er]div T (h n V T Ud) on duj, 
on dfl. 



(6) 



Au^ = in CI \ u) and in w, 

[u' n ] = h n — d n u~ on duj, 

[ad n u'„\ = [er]div T (h n W T u n ) on dui, 

du' = on dfl. 



(7) 



The main result of this work concerns the second order derivative. It is given is the following theorem. 



Theorem 2 Let fl be an open smooth subset of M d (d = 2 or 3) and let lo be an element of fld a 
with a C < a boundary. Let hi and h 2 be two deformation fields in Tt. Then the state Ud has a second 
order shape derivative u" d £ H 1 (fi \ZJ) UH'fw) that solves 

Au'd = in ft \ To and in to, 

[ u d] = {hi,nh2, n H - hi T .(Z?nh 2T )) [d n ud] - (hi, n \d n (ud)' 2 ] + /i2,n[dn(«d)'i]) 
+ (hi r .V/i 2 ,„ + h 2T .V/ii in ) [d n Ud] on dio, 

[trd n u'£\ = div T (h 2 , n [o-V T (««j)i] + h hn [aV T (u d y 2 ] + h lr .(£>nh 2T )[crV T u d ]) (8) 

- div r f(hi T .V r /i 2 ,„ + V T /ii jn .h 2T ) [o-VrUd]) 
+ div r [h 2n hi n (2D\i — HL) [a\7 T Ud]) on dio, 
u'J, = 0on9O. 
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Here, (ud)'i denotes the first order derivative of u in the direction of hi as given in ©, Dn stands for 
the second fundamental form of the manifold dm and H stands for the mean curvature of duj. The 
twin result concerning u n is an easy adaption of Theorem [2j Once the differentiability of the state 
function has been established, one can consider the objectives. In [I], we have shown the first order 
result. 

Theorem 3 Let £1 be an open smooth subset of M. (d — 2 or 3) and let u> be an element of ild 
with a C ' a boundary. Let hi and h.2 be two deformation fields in hi. The Kohn-Vogelius objective is 
differ entiable with respect to the shape and its derivative in the direction of a deformation field h is 
given by: 

CTl 

02 



DJ KV (uj)h = [o-} f — (\d n u+\ 2 - \8 n u+\ 2 ) + \V T u d \ 2 - \V T u r , 



(9) 



We now give the second-order derivative of the Kohn and Vogelius criterion. 



Theorem 4 Let Vl be an open smooth subset ofM. d (d — 2 or 3) and u be an element of fld with a 
C 4 ' a boundary. Let hi and h.2 be two deformation fields in hi. The Kohn-Vogelius objective is twice 
differ entiable with respect to the shape and its second derivative in the directions hi and \i2 is given 
by: 



D 2 J KV {iu){h 1 M) 



o-IVul 



(h lr .V T (fta,n) + h 2r .V r (/n,„) - h 2r .(Dnh lr )) 



On 



a\Vv\ 2 h hn h 2 , n + 2 aVv.(h hn Vv' 2 + h 2 , n Vv' 1 ) 



o~ ( 9 n (wn)i^2 + d n (u n )' 2 v' 1 — d n v[(ud) 2 - d 



(10) 



2 / v 



crd n (u n )". 



aid n v + 



where we have set v = Ud — u r< 



To investigate the properties of stability of this cost function, we are led to consider an admissible 
inclusion oj* to solve both (H}-© an d ©-(HI) in order to obtain the corresponding measurements /* 
and g*. It is obvious that the domain lu* realizes the absolute minimum of the criterion Jkv since, 
by construction, we can write Ud = u n in f2 and hence Jkv{^*) — 0- We will check that the Eulcr 
equation 

DJ KV (iu*){h.) = 0, 

holds. We will also prove that 



£ 2 JKv(^*)(h,h) = / <t\Vi 
Jn 



(11) 



Moreover, if h n ^ 0, then D 2 Jkv ( w *)(h, h) > holds. Nevertheless, (fTTj) does not mean that the 
minimization problem is well-posed. In fact, it is the following theorem that explains the instability 
of standard minimization algorithms. 

Theorem 5 Assume that to* is a critical shape of Jkv for which the additional condition u n = Ud 
holds. Then the Riesz operator corresponding to D 2 Jkv{^*) defined from H 1 / 2 (9u;*) with values in 
H^ 1 ^ 2 (duj*) is compact. Moreover, the minimization problem is severely ill-posed in the following 
sense: if the target domain is C°° and if X n denotes the n th eigenvalue of D 2 Jkv{^>*), then \ n = 
o(n~ s ) for all s > 0. 
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Theorem [5] has two main consequences. First, the shape Hessian at the global minimizer is 
not coercive. This means that this minimizer may not be a local strict minimum of the criterion. 
Moreover, the criterion provides no control of the distance between the parameter uj and the target 
u>*. The second consequence concerns any numerical scheme used to obtain this optimal domain w* . 
One has to face this difficulty and this explains why frozen Newton or Levenberg-Marquard schemes 
have been used to solve numerically this problem pj [T] . 

The paper is organized as follows. In a first section, we state some preliminary results. Some are 
well known facts in shape optimization and will be recalled without proof for the sake of readability. 
Some of them (e.g the derivatives of a Laplace-Beltrami operator and the tangential regularity of the 
solution to (|T|)- dZ|) along the discontinuity of the conductivity distribution) are less known and will 
be proved thanks to potential layer methods. Hence we will tackle the computations in Section [3] 
that we consider as the core of this work : it is essentially devoted to prove Theorem [5] After a first 
part where we prove the existence of a second order derivative for the state, we propose two distinct 
methods to find the boundary value problem solved by this second order derivative. The first method 
(subsection 13. 3p follows the lines of classical proofs of shape differentiability by differentiating the 
weak formulation of problem CP)-© and interpreting the result in terms of differential operator and 
boundary conditions. The alternative method (subsection 13. 4p consists in a direct differentiation of 
the boundary conditions. Finally, Section [4] is devoted to the analysis of the criterion, we establish 
Theorem [4] and Theorem [5] We will present their consequences on the stability of critical shapes. 



2 Preliminary results. 

2.1 Elements of shape calculus 

Before entering the proof of Theorem [2l we recall without proof some basic facts from shape opti- 
mization (see [6] for references). Let h be a deformation field in C 2 (Sl,R d ) with ||h||ca < 1. We set 
T t (h, .) = Id + th and denote by fit the transported domain Sl t = 2t(Q). To avoid heavy notations, 
we will misuse the notation T t instead of T t {h, .). 

Material and shape derivatives. Classically, in mechanics of continuous media, the material 
derivative is defined as being a positive limit. In our context, for any vector field h G 7i, we define 
the material derivative of the domain functional y = y(Q) at £1 in an admissible direction h as the 
limit 

y(n t )oT t -y(Q) 
y ( ' ^ = 1™ 1 ' ^ ^ 

Similarly, one can define the material derivative y(dQ, h) for any domain functional y = y(dQ) which 
depends on dil. Another kind of derivative occurs : it is called the shape derivative of ?/(fi, h). It is 
viewed as a first local variation. Its definition is given by the following 

Definition 1 The shape derivative y' = y'(f2;h) of a functional y(Ct) at ft in the direction of a 
vector field h is given by 

y' = y - h.Vy. (13) 

For more details on these derivations, the reader can consult 15]. 

Elements of tangential derivatives. We will need in the sequel to manipulate the tangential 
differential operators on a manifold. For the reader's convenience, we recall from some definitions 
and also some useful rules of calculus. 

Definition 2 The tangential divergence of a vector field V £ C 1 (M d ,M d ) is given by 

div T (V) = div (V) - DV.n.n, (14) 
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where the notation _DV denotes the Jacobian matrix of V. When the vector V £ C 1 (dfl,M. d ) is 
defined on dQ, then the following notation is used to define the tangential divergence 

div r (V) = div (V) - (DV.n).n, (15) 

where V stands for an arbitrary C 1 extension o/V on an open neighborhood of dCl. 

We introduce now, the notion of tangential gradient V r of any smooth scalar function / in 

Definition 3 Let an element f £ C 1 (9f2,R d ) be given and let f be an extension of f in the sense 
that f £ C 1 (?7) and f\an = f and where U is an open neighborhood of dft. Then the following 
notation is used to defined the tangential gradient 

V r / = V/| 9 n - V/.n n on dCl. (16) 

The details for the existence of such an extension can be found in [4]. Let us remark that these 
definitions do not depend on the choice of the extension. Furthermore, one can show the important 
relation 

/ V T /.P = - / / div T (F) , (17) 
Jan Jan 

for all elements / £ C l (d£l) and all vector fields F £ C 1 (5fi,M d ) satisfying F n = (F, n) = 0. 

Integration by parts on dil. In general, the condition above F n = is not always satisfied. 
We are then led to find another formula to extend the formula in the general case. The extension of 
this integration by parts formula to fields with a normal vector component involves curvature. 

First, we point out that the curvature is connected to the normal vector via the tangential 
divergence operator. Recall that the mean curvature of <9f2 is defined as H = div T (n). Making use 
of the form of div r (n) on the boundary, one shows straightforwardly the following statement. 

Proposition 1 Let Q be an open subset ofM 3 with a C 2 boundary. For any unitary extension Af of 
n on a neighborhood of dVt, one has 

div (TV) =H ondVt. 

Assume that the manifold dfl has no borders. IfF £ H 2 (dQ) 3 and f £ H 2 (dQ), then we have 

f V/.F + /div T (F)= / (V/.n + iJ/)F.n. (18) 
Jan Jan 

We assume now that the domain SI has a C 3 boundary. The simplest second-order derivative is the 
Laplace Beltrami operator; it is defined as follows (see OSHS]) thanks to the following usual chain 
rule. 

Definition 4 Let f £ H 2 (dfl). The Laplace- Beltrami A T of f is defined as follows 

A T f = div T (V T /) . (19) 

There is a relation connecting the Laplace operator and the Laplace-Beltrami operator. Let us denote 
by d 2 n f = (D 2 f. n).n where D 2 f stands for the Hessian of /. 

Proposition 2 Let 11 be a domain with a boundary dQ, of class C 3 . For all functions f £ H 3 (SY), it 
holds 

Af = A T f + Hd n f + d 2 m f 1 on Oil. (20) 
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We need to compute shape and material derivative of special vector fields: the outer unit normal 
vector n, the tangential gradient and the Laplace-Beltrami operator applied to a function. While the 
derivative of the normal vector is obtained by a straightforward calculus, we have to transport from 
<9f2t to dfl the Laplacc-Bcltrami operator and the tangential gradient in order to compute the other 
derivatives. 

Derivatives of the normal vector. We describe the material and shape derivatives of the normal 
vector. We will denote by n the gradient of the signed distance to d£l. This is an unitary extension 
of the unitary normal vector n at d£l which is smooth in the vicinity of dQ. This extension furnishes 
a symmetric Jacobian Dn that satisfies Dn n = on d£l. The direction h will be supposed to be in 
C 2 (M d ,R d ) or in C 2 (dfl,R d ). 

Proposition 3 The material derivative n of the normal vector n at in the direction of a vector 
field h G C 1 (M d ,M d ) is given by 

ri = -V T (h.n) + Dnh r , 

where h T = h h.n n. 

Concerning its shape derivative defined as n' = (dtiit)\t=o where n t is any smooth unitary extension 
of n to dflt, we obtain. 

Proposition 4 The shape boundary n' in the direction of h is given by 

ri = -V T (h.n). 

Derivative of the tangential gradient. For / G H 3 (9f2), we compute the material derivative of 
V T /. We first compute the difference V T f — V/. 

Proposition 5 For all functions f G C 2 (M 3 ) and directions h G C 2 (cMl,R 3 ), one has 



V T / = V/ + (D 2 fh) T - V/.n h - V/.h n 
Proof of Proposition [5} We differentiate V/ and V/.n n and obtain 

Vf = Vf' + D 2 fh 

while 

V/.n n = V/.h n + V/.n h + V/'.n n + (D 2 /h).n n. 
The two former equations give the desired result. I 
Derivative of the Laplace-Beltrami operator. Now, we want to compute the material derivative 



A r /. We begin to study how to transport the Laplace-Beltrami operator when one works on dQf 
Let A Tjf denote the Laplace-Beltrami operator on the manifold dilf To compute the derivative of a 
Laplace-Beltrami operator, we need the following proposition that we quote from [20] . 

Proposition 6 Let f G R 5/2 (R d ), then 



on 



(A T , t /)oT ( 7r (t) 



Oil 



C{t) V(/oT t )-(S(t) n).V(/oT t ; 



(21) 
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In the former proposition, we set 

T(*) 

7r(*) 

B(t) 
C(t) 

A straightforward computation gives 

7 '(0) = 

y T (o) = 

fi'(O) = 
C"(0) = 



det DT t , 

7 (t)||( J DT t - 1 ) T .n|| R ,, 

7T (0^(rr 1 )( J D(T t )- 1 ) T . 



div r (h) , 

div T (h) = div T (h T ) + Hh n , 
2(Dhn).n/ - (Dh + (Dh) T ), 
div r (h) / — (Dh + (Dh) T ). 



(22) 



(23) 



Theorem 6 Let f £ P(K d ). T/ie material derivative of A T f in the direction h is given by 



A T f = A r /+V r /.V r [div r (h r )] + V T (Hh n ).V T f - div T Dh+ (Dh) T j V T / 



(24) 



Proof of Theorem [6] : Formula (|24f is shown in a weak sense. For each test function <f> e C°°(dft), 
there exists an extension <f> G V(R d ) such that d n <p = 0; this can be done by extending ^ as a constant 
along the orbits of the gradient of the signed distance function to dft and the use of a cut-off function. 
For / G D(M d ), we set 

f (A T , t /)oT t -A T / 

A(i) = / 7 T (t) 0. 

Jan r 



After an integration by parts on <9f2, we obtain: 

1 - 7r(*) 

/an * 



A(t) = 



/ 1 J tW (M)or, <f,+ f ^ f(A T)t /)oT t + iv T /.V T( 
Jan 1 ian 1 V 1 



l~7r(*) 



an 



(A T ,t/) ° H 



/ 7 ( t v -/ ~ C W v (/ ° T t)} ■ W + [(B(f)n.V(/ o T t )] C(i) n. V<£) . 
Jon 1 K ' 



Since d n </> = and C(0) = /, we get 

/" l-7r(t) 

/an * 



A(t) = 



(Ar, t /)oT t 



V T (/-/oT t ) 



v T <; 



an 



C(0) - C(t) 



an 



V(/oT t ).V r( 



When t — > 0, it then comes 



/ = - / i T (t)A T f<f> + V r /.V r + (c"(0).V/) .V T <£, 

Jan Jan v ' 



an 



an 



A T f - div T (h) A T /j cf> + \Dh + (Dh) T - div T (h) I) Vf.V T <t>, 

A T f - div r (h) A T f + div T (div r (h) V T /) - div r [ (7-Dh + (£>h) T ) V/ 



Expanding the double divergence term, we obtain: 



A T f = A T f + V r /.V T div T (h) - div r \({ph+ (I?h) T ) v/) ^ . 
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In order to explicit these derivatives, we let appear the curvatures of <9f2 by means of 

V r /.V T div T (h) = V r /.V T [div r (h T ) + Hh n ] , 
and this ends the proof of the theorem ([24]) . 



3 Existence of the second order derivative of the state. Proof 
of Theorem [2l 



The section is devoted to prove Theorem O We follow the usual strategy to derive existence in shape 
optimization. In section [221 we will write the weak formulation of the problem, then transport it on 
the reference domain, pass to the limit and obtain existence of the material derivative. In a second 
time, we will seek a boundary value problem solved by the material derivative. This will provide 
a characterization of the second order shape derivative. Two strategies, that we will detail, are 
possible: the first one explored in section [3731 consists in working on the variational formulation while 
the second one uses the tangential differential calculus by differentiating the boundary conditions. 
This last approach will be presented in section l3~4l The computations that will be made in subsections 
13.31 and 1 3.41 require some regularity of the traces of the state Ud on the interface of discontinuity duj. 
For the sake of readability, we postponed in subsection 13.51 all the needed justifications. 

3.1 Preliminary results. 

In the sequel, we will use some technical formulae. To preserve the readability of the proof of the 
main result, we state them in this paragraph. The tools needed for proving these results can be found 
in [2D]. Given a smooth vector field h, we denote 



A h = Dh + Dh T - div (h) / 



We begin with the following formula. 



Lemma 1 It holds: 



\7u.A h \7v = V(h.Vu).Vv + V(h.Vv)Vw - div ((Vu.Vw)h) 



(25) 



Given two smooth vector fields hi and h 2 , we set 



21 = Dh 2 A hl + A hl Dh 2 T - A hl dW (h 2 ) - {A hl )'(h 2 ), 



(26) 



and 



b = (h 2 .X7u)A hl X7v + (h 2 .Vv)A hl Vu - {{A hl Vu).Vv)h 2 . 
Here, the notation (Ahi)'(h 2 ) stands for the matrix defined by its elements 

((^ hl )'(h 2 )) M =V(((A hl )')M)-h 2 



Lemma 2 One has: 



Vm.SIVw = div (b) - (h 2 .Vu)div ((A hl Wv)) - (h 2 .Vv)div ((A hl Vu)) . 



(27) 



We need the following crucial result 
Lemma 3 If u is harmonic then 



div (A hl Vu) = A(hi.Vu). 



(28) 
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Proof of Lemma [3] For any harmonic function u in f2 and for every test function <fi € T>(£1), we 
can write 

/ VuV</> = / A h Vu\7cf) 
Jn Jn 

then 



/ Aii (j) = / div (A h Vw) 



Since w = u' + h.Vu and since u' is harmonic in Q, we obtain the result. I 
3.2 Proof of existence of the second order derivative. 

We follow Hettlich and Rundell [8] and Simon [19] to define the second order derivative of an op- 
erator with respect to a domain. We compute the second derivative by considering two admissible 
deformations hi, h 2 £ H that will describe the small variations of doj. Simon shows that the second 
derivative F"(duj; hi, h 2 ) of F{duo) is defined as a bounded bilinear operator satisfying 

F"(5w;hi,h 2 ) = (>(&,;; hi)) 'h 2 -F'(&^hih 2 ) (29) 

For more details, the reader can consult the appendix in page 613 of [8]. 

Let us begin the proof. Let hi,h 2 g Ti be two vector fields. The direction hi being fixed, we 
consider ui^ 2 the variation of u\ with respect to the direction h 2 . We recall from [T] that the material 
derivative it'i of u in the direction hi satisfies 



Jn Jn 



Let </> 2 : Q i— > 17 be the diffeomorphism defined by 2 (x) = x + h 2 (x) and we set -02 = 4> 2 1 ■ Setting 
^h 2 = { x + h 2 (x), x e w}, f2h 2 = {x + h 2 (x), x S 0} = and Oh 2 = cr ° 02, we get 

o-h 2 Vui jh2 -Vw = / ffh 2 Vw h ,.4 1 Vi' (30) 
Jn h2 

where i*h 2 is the solution of the original problem with Wh 2 instead of w. Making the change of 
variables x — <f> 2 (X), we get the integral identity on the fixed domain il : 

J crVZ 1M .^DMDip2) T det(D(j) 2 )^Vv = J aVu h2 . (Dtp 2 A^(D^ 2 ) T detpfc)) V« (31) 

with the notations u = u o <p 2 and Ah t = o <f> 2 . Since the material derivative ii\ of u with respect 
to the direction hi satisfies 

/ <tVu\.Vv= j <7\7u.A hl Vv, 
Jn Jn 



the difference of (|30j) and (131]) gives 

J crv(ii 1M -u^j.Vv = J aVZ 1M . (l-D^ 2 (Dip2) T det(D<p 2 fj Vw 

+ J aVu h2 . [Dip2A\ 1 {D^ 2 ) T det(D(f> 2 ) - A hl ) Vv + ^ (Vu h2 - Vu).A hl Vv. 
We quote from [13] and [8] the following asymptotic formulae 



div^oiU = odlhiWl*), 



IC 2 

WD^D^f" fet(D<t>i) - 1 + A K \\ X = Odlhillg,), 



ll-D^^hi (Dtp 2 ) T det(D0 2 ) - A hl + Z)h 2J 4 hl + A hl (Dh 2 ) T — div (h 2 ) A hl - (AJ'MU = O(||h 2 ||c0- 
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Making the adequate substitutions, we easily check that the material derivative of u\ with respect 
to h 2 exists. This derivative, denoted by u\, satisfies 



/ aVui.Vv dx= I (7 [Vui.Ah a Vw + Vii 2 .A hl Vv - VuMVv] 
Jn Jn 



where 21 is defined in ([26 



(32) 



3.3 Derivation of (JSJ) from the weak formulation. 

We want to make explicit the problem solved by (u')' . To achieve this, we should write the right 
hand side 

F = I a[Viii.A h2 Vv + Vu 2 .A bx Vv - VuMVv] , 

as the sum of an integral with Vt> in factor and an integral of a divergence to identify the jump 
conditions on du). To that end, we will use algebraic identities that involve second order derivatives 
of u, in and of the test function v £ 2?(f2). Using Lemma [TJ we obtain: 



aVui.Ah 2 Vv = I a 
n ~ Jn 



erV-U2-^4hi Vv = / a 
n Jn 



V(h 2 .Vui).Vu + V(h 2 .Vw)Vui - div ((Vui.Vw)h 2 ) 



V(hi. Vu 2 ). Vv + V(hi.Vu) Vu 2 - div ((Vii 2 .Vw)hi) 
Concerning the remaining terms, we use Lemma [5] to get 

/ aVuMVv = / a div ((h 2 .Vu)A hl Vv + (h 2 .Vv)A hl Vu - (A hl Vu.Vv)h 2 
Jn Jn 

- a (h 2 .Vu)div (A hl Vv) + (h 2 .Vu)div (A hl Vu) 
We apply Lemma [3] and gather the expressions obtained for F. 



f = I a V (hi.Vu 2 + h 2 .Viii) .Vu + V(h 2 .Vw).Vui + V(hi.Vv).Vii 2 

f / a div ((A hl Vu.Vv - Viii.Vv)h 2 - (Vu 2 .Vw)hi) 
Jn 

(h 2 .Vu)A(hi.Vu) - div ((h 2 . Vv)A hl Vu) - V(h 2 .Vu).A hl Vv 

Using (|23|) . we remove the dependency on Ah 1 Vv: 

V(h 2 .Vu).A hl Vv = V(hi.V(h 2 .V«)).V« + V(hi.W)V(h 2 .Vit) - div ((V(h 2 .Vu).Vu)hi) 
Therefore, we write F — F\ + F 2 where 



(33) 



F 2 



.Vv, 



V (hi.Vu 2 + h 2 .Vwi) - V(hi.V(h 2 .Vu)) 
V(hi.V«).V(u 2 - h 2 .Vu) + V(h 2 .V«).Vui + (h 2 .Vv)A(hi.Vu) 
J a div ({A hl Vu.Vv - V«i.V«)h 2 + (V(h 2 .Vu).Vu - Vu 2 .Vv) hi - (h 2 .Vv)A hl Vu 



(34) 
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The connection between second order material and shape derivatives is given by: 

ill = «)' 2 + hi.Vn 2 + h 2 .VMi - hi.V(h 2 .Vu), 
incorporating this expression in (|34p , we rewrite (|32[) as: 



VvEH^fl), [ (tV(u' 1 )' 2 X7v = F 2 . 
Jn 



(35) 



Testing it against v G T>(Q, \ dui), we get A(w' 1 )2 = in £1 \ ZJ and in u>. We now deduce the jump 
conditions for (u'i) 2 . To obtain the jump of the potential, we simply write that u\ € H^fJ), hence 
[Hi] = on duj and then 

[K) 2 ] = -hi.V«2 ~ h 2 .Vui. 



To express the jump of the flux, we then apply the Gauss formula in (|35p to get 

- / [crdn(u'i)' 2 ]v = F 2 . 

J duj 



(36) 



The second term F 2 contains all the jumps of the flux on the interface du. 

A simplified expression of F 2 . To get a simplified formula for F 2 under a boundary integral, 
some lengthy but straightforward calculations are needed. We summarize the result by means of the 
following lemma 

Lemma 4 One has: 

F 2 = div T (2h 2n hi n Dxi [aV T u] — h 2n n.Vhi^ n [<tV t u] + h 2n h lT .Dn n [ctV t m]) v 
+ / div r (h lT .V T (h 2 ,n) [ctV t m] - hi, n h 2 , n H [erV T u]) v 

J du) 



(37) 



div r /i 2j „ aVrUi 



div T hi 



Proof of lemma First, write : 

/ aV(hi.Vv).V(u 2 - h 2 .Vu) = a x / V(hi.Vw).Vu 2 + a 2 / V(hi.Vu).Vu 2 
Jn Jn\u ~ Ju> 



[ad n u' 2 ](hi.Vv) 

Note that the normal vector is oriented from lj to £1 \ ZJ. In the same spirit, we write 

V(h 2 .Vw).Vtti + (h 2 .Vu)A(hi.Vu) = V(h 2 .Vw).V(ui - h x .Vu) + div ((h 2 .Vw).V(hi.Vw)) . 
By a argument of symmetry, we then can write: 

/ crV(h 2 .Vi>).V(?ii - hi.Vu) = - / [ad n u[](h 2 .Wv). 

J n J duj 

To drop the dependency in , we use (I25[) and get after expansion: 
div ((A hl Vu.Vu) h 2 ) = div ((V(hi.Vw).Vu + V(hi.Vu)Vw) h 2 ) - div (div ((Vu.Vw)hi) h 2 ) ; 
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div ((h 2 .Vv)A hl Vu) = S7{h 2 .Vv).A hl Vu + (h 2 .Vw)div (A hl Vu) 

= V(h 1 .V(h 2 .Vw)).Vu + V(hi.Vu)V(h 2 .Vu) + (h 2 .Vu)A(h 1 .Vu) 

-div((V(h 2 .Vu).Vw) 
= V (hi.V(h 2 .Vi;)) .Vu + div ((h 2 .Vi;)V(hi.Vu) - (V (h 2 .V v) .V u) hi) . 
After integrating by parts, we conclude thanks to the state equation and obtain 

f ffV(hi.V(h 2 .Vv)) .Vu = - [ (hi.V(h 2 .Vu))div(crVit) =0 
Jn Jn 

We substitute the shape derivative u' to the material one u\ 



F 2 = - f [cr9 n ui](h 2 .V«) + [o-d n ?4](hi.Vu) - f a div (div ((Vu.Vu)hi) h : 

+ J a div H (V(hx.Vu).Vu)h 2 + (V(h 2 .Vw).Vu) hi) - ((Vi4.Vu)hi + (V« x .Vi;)h; 

First, we use the continuity of the flux on du>, then we integrate by parts on du and finally we 
incorporate the expressions of the jumps of the shape derivatives u' to obtain 

/ crdiv (hi. (V(h 2 .V«).Vtt)) = - / [aVu.V(h 2 .Vv)] hi, n = - / [aV T u]h^ n V T (h 2 .Vv) 

Jn ^ ' Jdu Jdu 



div T ([o-V t m] hi :n ) h 2 .Vw 



du 



h 2 .Vv. 



aVu[.Vv 


= ~ I h 2:7l 


ad n u[ 


d n v - / h 2y7l 


aV T u[ 


V T v 




J du 




J du 







d n v 



du 



div r /i 2j „ crV T it'i 



This leads to a simplified expression for F 2 : 

F 2 = - J a div ^div ((Vu.Vu)hi) h 2 + ((Vui.Vv)h 2 + (Vu 2 .Vu)h] 

Let us study each term of this sum. Using Gauss formula and integrating by parts on the manifold 
du , we obtain 

/ a div (Vui.Vu)h 2 ) = - / h 2 , n 

L V ' Jdu 

= - I h 2 . n a9„u'i 

J du 

By symmetry, we also get: 

a div (vu' 2 .Vv)h 

We now turn to the term with a double divergence. We first write it as a boundary integral thanks 
to Gauss formula as 

/ a div (div ((Vu.Vu)hi) h 2 ) = J ft, 2: „div (hi [a(Vu.Vv)] ) , 
then, we use (fill) to introduce the tangential operators 
/ cr div(div((Vu.Vw)hi)h 2 ) = / h 2 . n div T (hi [cr(Vu.Vw)] ) + / h 2 . n D{h.i [a(Vu.Vvj\ )n.n. 

Jui ^ 1 Jdu ^ ' Jdu 
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) = - / hi, n 


ad n u' 2 


d n v 4 


/ div T hi^ n 


aV T u' 2 


)• 


1 Jdu 






Jdu \ 







We study each of these terms. We start with the one involving tangential derivatives: we expand the 
tangential divergence to incorporate the jump relation for the state u. 

div r (hi [er(Vu.Vt>)]J = div r (hi) [a(Vu.Vv)] +hi.V T [aVu.Vv] 

= div T (hi) [crV T u] .V T v + hi.V T [aS/ T u.V T v] . 

Then, the first term becomes: 

/ ft.2,ndiv T ( hi [cr(Vu. Vv)] ) = / h 2tn div T (hi) [crV T u] .V T v + / ft. 2jn hi.V T [a\7 T uV T v] . 



We use the integration by parts formula (fTB)) to get: 



ft<2,ndiv T ( hi [<t(Vm.Vi;)] 
= / hi :n Ii2 : nH [ctV t w] .V t v — div T (div T (hi) hi, n [cV T u]) v — div r (hili2, n ) [crV T u] .V T v 
= J div r f (div r (hi/i 2 ,„) - div r (h x ) h 2 , n - hi^ n h 2 ^H^ [erV r u]^ v. 
Expanding 



div T (div r (hi/i 2 ,„) [<rV T u]j v = div r (div r (hi) h 2 , n [crV T u] + hi.V T (h 2 , n ) [<tV t u]) i; 

= div T (div r (hi) h 2 . n [crV T u]) v + div r (hi T V r /i 2 ,„ [erV r it]) v, 

we obtain the new expression: 

J /i2,ndiv T (hi [er(Vu.Vu)]) = J div r ((hi T V T /i2,n - hi, n h 2 , n H) [crV T u]j v. 

Now, we consider the term involving normal components. We have 

n.D(hi [crV(Vu.Vu)] )n = n.V(/ii,„ [<rVii.Vu]) - [ctVu.Vw] hi T .L>nn 

= n.V(/ii, n ) [crV r u] V T w + /ii i „n.V([crVu.Vw]). 

Then, we get 

/ h 2 , n D (hi [cr(Vu.Vw)] )n.n = / /i2,n n -V(/ii, n ) [crV T u] .V T u + h 2tn hi, n n.V([crVu.Vv]) 
= / -div r (/i2,nn.V(/ii, n ) [£7V T u]) v + h 2>n hi !n ii.'V([a'Vu.Vv]). 

J duj 

A straightforward calculus leads to 



(38) 



(39) 



n.V([crVu.Vv]) = n. 



aD z uVv 



■ D 2 V [crVu] 



n. yd n v 
d n v 



<rD 2 u 



oD 2 u 



V T v + D 2 v [a\7 T u] 



d 2 u 








+ n. 


crD 2 u 


dn 2 







D 2 v [aV T 
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where D 2 u is the Hessian matrix of u. From (|20[) and from the jump conditions for the state u, we 
deduce that 

d 2 u 

= " [<TAtU] ■ 

When one differentiates the relation expressing the continuity of the flux for the state along the 
tangential direction V T u, one gets ([B], p 235): 

= V[ad n u].V T v = [aD 2 u]V T v.n + [aVu\.(DnV T v). 

In the same spirit, it comes that 

V<9 n v.[crV T u] = D 2 v[aV T u].n + V«.(Dn[crV r u]). (40) 

Since Dn is a symmetric matrix and Dnn = 0, one checks Vv.(Dn[aV T u]) = [erVu].(Dn V T v). 
Then 

n.V([crVu.Vw])) = - [crA T u] d n v - 2Dn [aV T u] .V T v + [aV T u] V T d n v 
We integrate this expression on du; and obtain after some integration by parts: 

/ /i2, n hi,nn.V([crVu.Vu]) 

J du 

= - h2,nhi,n [crA T u] d n v + / h 2 . n h 1:n [<jV t u] V T d n v -2 h 2 , n hi, n Dn [aV r u] .V r t> , 

J du J du J du 



du 



h2, n hi, n [crA T u] + div r (h2, n h>i,n [fxV r w]) d n v + 2 / div r (/i 2)n /ii )7l Dn [aV r u]) u. 



Hence 



/i 2l n^(hi [o-(Vu.Vu)] )n.n = — ^2,n^i,n [o"A r u] + div r {h 2i nhi. n [<rV T u]) 



du 



+ / div T (2ft,2,n^l,n£'n [cV T u] — /l2,n*l.V T /li jrl [aV T «]) w. 

crdiv ( div ((Vu.Vw)hi) h2 ) = / h 2 , n hi >n [crA T n] + div T (h 2tn h\^ n [cV T n]) 3 n t> 

— / div r (2/i2 i „/ii i „Dn [crV T u] — h 2 ,n n -^hi >n [crV T u]) u 

- / div r (hi r .V r (/i 2 ,n) [crV r u] - hi jn h 2 , n H [o\J T u]) v. 
Gathering all the terms, we write F 2 as: 

F 2 = J div r ^2/i 2 , n /ii, n -Dn [erV r u] + (hi T .V T (/i 2 ,n) - h 2 , n n.Vhi^ n - h\ tn h 2 , n H) [aV T u]J v 
^div r (^h 2 , n aV T u[ ^ + div r (^h\, n (jV t u' 2 

- I (h2,nhi, n [aA T u] + div T (h 2>n hi in [aV T u]) J d n t; 

- / (/ii,ndiv T (/i 2 , n [cVtu]) + h 2 , n div T (hi jn [crV T u]) ) d n v. 
We end the proof after expanding the tangential divergence of the last term of F 2 . 
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Let us return to the weak formulation (|36[) of the derivative. By identification, we get 



[cr9n (1*1)2] = div r ^/i2,n cr'V T u' 1 ) + div T yhi^n aV T u 2 J — div T (/i2,n^i,n(2£ ) n — H I) [o"V T ii]) 
- div r (h lT .y T (h 2 ,n) [crV r u] - /i 2 ,„n.V/ii,« [erV T u] + fc. 2) „h lT ..Dnn [aV T u]) . 
It remains to compute the jump of the flux for the second order derivative. Since 

<2 = («i)a-«Dh 1 h a (41) 

where u' Dhi ha is the first shape derivative of u in the direction of the vector field Dhi h 2 . Thanks 
to ([6]), we can write the jump under the form 



[<7<9 n u" 2 ] = [cr5„(ui) 2 ] - Wd n u' Dhih2 ] = [cr9 n ("i) 2 ] - div T ( J Dhih2.n[crV T u]) 



(42) 



Let us split the field h 2 in two parts: Dhx h 2 .n = ft, 2 „n..Dhi n + Dhx h 2r .n. In the spirit of (j4H)l . 
we obtain 

Z3hih 2T .n = V T /ii, n .h2 T — hi T .Z?nh 2T . (43) 
Thanks to (139j) . the jump [adnu' Dili h ] then can be written under the form 

[cr<9„u^ hlh2 ] = div r ((/i 2 ,„n.V/ii,„ + V T ft,i irl .h 2T - hi T .L»nh 2T )[crV T u]) . 
Gathering all the terms, simplifications occur and we get: 



[ad n u" 2 ] =div T I Ji2,n cVt-Uj 



hi 



aV T u' 2 



div T ((hi T .V T /i 2 ,n + V T /ii,„.h 2r ) [erV T u]) 



- div T (h 2 , n hi^{2Dn - HI) [crV T u]) + div T (hi T .L»nh 2r )[aV T it]) . 
To get the jumps of the potential, we use l|4ip and obtain 



1,2 



2 



K) 

"/l2,n 



Dhi h 2 
-fa, 



= -hi. 



<9 n it 2 



-h 2 .[Vui]- 



hi T . 



V T u 2 



/i 2 ,„n. [V(hi.Vu)] - fti,„n. [V(h 2 .Vu)] 



l Dh 1 h 2 



Thanks to the jump of the potential for the first order shape derivative given in ([6]), it comes that 



V r iti = -h 2r . [V(hi.Vu)] and h iT . V T u 2 = -h lr . [V(h 2 .Vu)] 



and then: 



d n u[ 



d n u'n 



/i 2 ,„n. [V(hi.Vu)] +h lr . [V(h 2 .Vu)] 



L Dh! h 2 



(44) 



Computing the other jumps that appeared in the former expression, we get 



[V(h 2 .Vu)] = (Dh 2 ) T [Vu] 



D 2 u 



D 2 u 



h lr . [V(ha.Vtt)] - n.L>h 2 h lr [d n u] + h 2 . n h lr . 
h 2 . n n. [V(hi.Vu)] = /i 2 ,„ [dau] n.Dhin + /i 2l „/ii, n n. 



n + hi T . 



D 2 u 



n + hi .„n. 



h 2r . 

D 2 m 



hi r . 



*-Dhih 2 



= -Dhi h 2 .n [<9 n u] = - (/i 2jn n.Dhin + n.Z?hih 2r ) [<9 n u] 
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With the help of formula (|4"5)) . we obtain: 
-h 2 , n n. [V(hi.Vu)] +h lr . [V(h 2 .Vi 



l Dh 1 h 2 



(V T /ii, n .h 2r + V T h, 2 ,n-hi r ) [5„u] 



2hi T .Dnh 2T [<9 n u] + hi T . 



D 2 u 



D 2 u 



D 2 u 



n = [A r w] - i? [<9 n it] = -ff [9 n u] 



hi r . 



D 2 u 



h.2 T = hi r .D([VM])h2 r = hi T .D([d n u] n)h2 T = [d n u] hi T ._Dnh2 r . 



Finally, we gather the results of these computations to write 



d n u[ 



hi. 



d n u' 2 



(V T /ii,„.h 2T + V T /i 2 ,n-hi r ) [d a u] 
{h2, n hi, n H - hi r .Dnh 2T ) [d n u]) 



(45) 



3.4 How to recover (JH]) by formal differentiation of the boundary condi- 
tions. 

The aim of this section is to retrieve the expression of the flux jump [ad n u"] by computing the normal 



derivatives of each of the expressions [crVu'].n and div r (/ii,„[crV r u]) . Since 

[crVu'J.n = div r (h l n [aV r it]J — hi jn [crA T u] + V T /ii, n .[erV T M], 

then, we get 

[cjVit'].n = hi^ n [aA T u] + hi^ n [aA T u] + V T /ii,n-[c r V r u] + V T hi. n .[aV ' T u]. (46) 

In order to avoid lengthy computations, we shall concentrate on each normal derivative appearing in 
the above formula. Some of the results are straightforward and their proof will be left to the reader. 
Combining propositions © and (O, we conclude that 



V r /ii,„ = - V T (hi.V T /i 2 ,n) + (£ ,2 /ii,„.h 2 ) T - V/ii,„.h n - V/ii,„.n h. 
In the same manner, we also get 

[ctVi-m] = [<tV t m 2 ] + ([c-D 2 «]-h 2 ) T — [erV r u].ri n — [erV T u]n h. 
Hence, we can write 



hi.n = h 2 .Vh „ - V T /i 2 ,„.h lT . 
It remains to simplify the terms A = (Z? 2 u.h 2 ) T .V T /ii in and B = [crV T u].(Z? 2 /ii.„.h 2 ) T . We obtain: 

A = -[aV T u].(Dn V T h^ n )h 2 ,n + [crA r u] V r /ii,„.h 2T , 

B = (D 2 /i2,„.h 2r ).[crV r u] + V r (9 n /ii i „).[crV T u]/i 2i „ - [ffV T «].(DnV T /ii, n )/i2,„. 

We tackle the computation of (d n u')' . For the sake of clearness, we subdivide the work in several 
steps. 

First step. We compute div T (/ii, n [crV r u]). We expand: 



div T (h lyn [aV T u\) = h ltn [aA T u] + V T hi in .[aV T u], 



hi yn [tjA T u] + hi, n [aA T u] + V T /ii !T! .[erV r u] + V T /ii ir! .[(7V T u]. 
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Hence, after substitution, one gets 

div T (hi^laVru]) = div T (hi^ n [<7V T u' 2 ] + (/i 2 ,, l 9„/ii i „ - V r /i 2i „.hi T )[erV T u] 

+ 2[aA T u]V T h 1 . n .h 2T - d n h 1 ^aV T u}.(Dnh 2T ) + [aV T u].(D 2 h 1 . n .h 2T 



2h 2 ,n[^ru].(DnV T h 1 ^ n ) + <jA t u - [aA T u 2 ] . (47) 



Second step. We compute [credit']]. From the expression of n, we get after some straightforward 
computations: 



[acUi] = [ad n (u[)' 2 } + ([aD 2 ^] h 2 ).n + [aV T u[].{Dn h 2r - V T /i 2 ,„). (48) 

Third step. We compute ad^Ui)^. From the jump condition on the flux of the derivative (JB]) and 
471) and (SHI), we obtain: 



[0-9(^)2] = div r (hi,nWV T u 2 ] + {h 2l nd n hi, n - V T h 2 . n .Yi lT ) [aV T u]j + 2V r /ii, n .h 2T [erA T u] 

-([aD 2 ^] h 2 ).n + [o-V T u[]. (V T ft, 2: „ - Dnh 2r ) - <9„/ii,„[crV T u].(D n h 2T ) 
+ (D 2 h 1 , n h 2r ).[aV T u] - h 2 . n (Dn [aV T u})£.\7 T h 1}7l . 

Taking account of the following calculation, 

-([cr£>V 1 ]h 2 ).n+ [aV T u' 1 ].V r h 2 ,n = - (h 2 . n {aD 2 u' 1 ]n+ [aD 2 u[] h 2r ) .n + [crV r u' 1 ].V r /i 2 , n, 

= h 2 , n ([aA T u[] + Hiadnu'^ + [au[}.V ' T h 2 , n - ([aD 2 u'l] h 2T ).n, 
= div r (h 2>n [aV T u[fj + H^.^ad^] - ([aD 2 u[] h 2T ).n; 

it comes 

[o-a n (ui) 2 ] = div T {hi. n [aV T u' 2 ] + /i 2 ,„[crV T u[] + {h 2tn dji^ n - V T h 2 ^ n .h lT ) [erV T u]) 
+2[crA T u]V r /ii,„.h 2T + Hh 2 , n [ad n u[] - {[aD 2 u[} h 2T ).n 

'[o-VrU^] + d n h ltn [o-V T u]) . {Dn h 2r ) + {D 2 h hn h 2r ) . [crV T u] 



-2/i 2 ,„V T /ii, n .(£>n [crV r u]) + h hn [[aA T u] - [aA T u' 2 ]j . (49) 

This formula remains hard to handle. To get a more convenient one, we decide to derive tangentially 
to the direction h 2 the boundary identity 

[ad n u[] = hi tn [aA T u] + V T h ltn .[crV T u\. 

This leads to: 

([aD 2 ^] h 2r ).n+(Dnh 2T ).[crV r u' 1 ] = V r /ii,„.h 2T [crA T u] + h 1 , n V T [aA T u}.h 2T 

2 (50) 
+ (D hi. n h 2r ).[crV T w] - <9 n /ii,„[o-V T u].(L>nh 2r ) + [crA r u]h 2T .V T /ii,„. 

From (|24p and subtracting (|5"0")) from (|4"9")> . we can write 

[0-<9 n (ui) 2 ] = div r (/ll.nlo-VrUy + /l 2 ,„ [<7 V r U^] + (/l 2 ,n<9 n ft.l : „ - V T ft, 2i „.hi T ) [crV r u] 

+div T (hx, n h 2 ^{Hl - 2L>n).[erV T u]) - /ii, n (V T [oA T u].h 2r + A r [crV r u].h 2T 



-/H, n V T div T (h 2r ) .[<rV T tt] - div T Dh 2 + (£h 2 ) J [ctV t m] 
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From (|24p . we obtain 



[aA T u] = [aA T u] + V r div r (h 2r ) .[<tV t u] + V T (Hh 2 ^ n ).[aV T u] 

(51) 



div T ^(jDh 2 + (Dh 2 ) T ) [aV T u] 
and using the relation between the material and shape derivative, we get 



[aA T u] = [aA T u'} + V ({aA T u}) .h 2 and [aA T u] = [aA T u'] + A T ([crVu].h 2 ) . 
Injecting these relations in (|51[) and applying them for h 2T , we get 

A T ([crV r u].h 2T ) + V r div T (h 2r ) .[a\7 T u] = V T [crA T w].h 2r + div r (j^(Dh 2 + (L>h 2 ) T ) [aV T u] 
This last fact allows us to conclude. 

3.5 Justification of the formal computations. 

We have to justify rigorously that the right-hand sides of ([5]) , , © make sense. They involve 
tangential derivatives of u n and Ud along the interface duo up to the order three. The existence of 
these derivatives is not clear a priori since the gradient of the solution has a discontinuity along this 
interface. Our first aim is to precise the tangential regularity along the interface dui of the solution 
u of ([T]) with either Dirichlet or Neumann boundary conditions. 

We should access to the trace of u on the interface dui. Any numerical discretization needs also 
to compute the state, its derivatives with respect to the shape and the normal derivatives along the 
interface du>. To that end, we introduce for any a G H^idw) and (3 G H~^ 2 {duj) the following 
boundary value problems 

Av — in CI \ ZJ and in lu, Av = in CI \ 57 and in ui, 

ID) { r = "° n f' and ^ r /! = r n f' (52) 

' [<JO n v\ = p on ouj, I [(70 n v\ = p on olu, 

v = ft on dCl. d a v = gi on dCl, 

where G H^ 2 (dQ) x H- 1 / 2 ^). Note that for a = 0, /3 = and (A, 31) = (/,p) then (u rf ) 

and u n solve respectively (D) and (N); furthermore the choice of 

[a] 

a = h n d n u and [3 = [a] div r (h n V T u) (53) 

02 

leads to ^ and ([7]) when we take 5) = (0, 0). 

Existence of solutions to (D) and (N). To study these problems, we use the integral rep- 
resentation in terms of layer potentials. In a first step, we recall some definitions. The Newtonian 
potential T is defined as: 



^-ln(|a; - y\) if n = 2, 

1 1 

- — 1 r if n = 3. 

47T \x — y\ 
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The integral equations applying to direct problem will be obtained from a study of the classical 
single- and double-layer potentials. We begin to introduce the following operators 



T(x,y)u(y) da{y); 



Sdujdnu(x) := / T(x,y)u(y) da{y); 



K d o,duju(x) 
K du)d Q,u{x) 



d n T(x,y)u(y) da(y) ; 
d n T(x,y)u(y) da(y) 



Note that all these operators have a smooth kernel since the boundaries dui and dfl are assumed to 
have no common point. We also denote 



Sn ■ u i-» Snu(x) : 
Kn : u i-» K n u(x) 



on 



K u u(x) 



T(x,y)u(y) da{y); 
d n r(x,y)u(y) da(y); 

T(x,y)u(y) da(y); 
d n T(x,y)u(y) da(y). 



on 



We now obtain some systems of integral equations to compute the state function and their shape 
derivatives. Since v is harmonic in Q \ 57 and for all x E dfl U duo, it has the classical boundary 
representation: 



1 



-v(x) = / d n T(x,y)v(y) 



on 



d n T(x,y)v(y) 



on 



T(x,y)d n v(y) 



Similarly since v harmonic in to, for all x G dcu we can write 

1 



-v (x) 



d n T(x,y)v{y) 



F{x,y)d n v(y). 



T(x,y)d Tl v(y). (54) 



(55) 



Ouj 



Let us denote by Vd the solution of the boundary values problem (D) in ([52]) . Let us show how 
to compute their restrictions and also their normal vector derivatives on the boundaries. Incorpo- 
rating the jump conditions, a straightforward computation leads to the following boundary integral 
equations 



1 

[l-Kdudn 



0"! 



02 + 01 
01 



Sdndu 



02 + 01 



-Sn 



1 



01 + 02 







a 


0i 








H 

01+02 


SdLudQ 









-02-KowOn 

where \i = [ct]/(<7i +02). Thanks to (|55p . the quantity (d n Vd) + is then given by 



(56) 



02 
01 



1. 



K„ 



1 

01 
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Concerning v n , the solution of the Neumann problem (N) in (|52p . the same kind of computations 
gives 



1 

n 1 + ^Ku. 



cr 2 + cti 



01 



a 2 + a i 



- 2 I + K Q 



0\ + 2 



(vn)\an 



5' 



(3 



(57) 



0"! 



01 + (72 



SdQdujgi 

Sngi 



Finally, the computation of (<9 n i>„)jt is given by 



S u (d n v n ) 



\du 



&2 I 1 . 
0"! 



Concerning the well-posedness of ([56")) , we can state the following result. 

Theorem 7 The linear system of integral equation (|56p /las an unique solution in ^^(dco) X 
H-V a (flfi). 

Proof of Theorem [7] Let A be the matricial operator defined on H 1 ' 2 (9o;) x H _1 ' 2 (9S1) as 



A = 



1 CTl 

-/ + fxK u ; SdCldu 

2 0~2 + CTl 



01 



01 + 02 



-<5n 



(58) 



The main argument of the proof is based on the Fredholm alternative. In a first step, we have to show 
that the adjoint operator A* is injective. Since the boundaries are bounded, the adjoint operator A* 
defined on R~^ 2 {duj) x H 1 / 2 ^) can be written under the form 



A* = 



-I + fiK* ^K* mduJ 



01 



02 + 01 



Sdi 



o-i 



udCl 



02 + 01 



-Sc. 



(59) 



Let (u,v) G R-^ 2 {duj) x H 1 / 2 (9S1) be in the kernel of A* . Consider the potential W defined for each 
x G K d by 

0i 



W(jb) 



02 + 01 \Jduj 



T(x,y)u(y) + / T(x,y)v(y) 



(60) 



In a first step, we show that W = 0. The function W satisfies AW = in R d \ (duj U 89) by 
construction. We check that W|an = from the equation corresponding to the second line of A* . 
By the properties of the single layer potential, [W] = on du>. Furthermore, it holds [adnW^] = on 
duo. Indeed, we can have ([TT]) 



d n W^ 



01 

01 + 02 



K* u + Kt DBui v 
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and 



hence, 



K* )v + K. 



dQdui v ' 



<Jid n W + - a 2 d n W- = a x I (-/ + ^K*)u + nK* mdu v J . 

This corresponds to the first line of A*(u,v). Then, W solves the Laplace equation |T]) with ho- 
mogeneous Dirichlet boundary conditions. By the uniqueness of the solution, we get W = in 

n. 

In a second step, we deduce that u — v = 0. Since W = in f2, we see that [d n W] = on dui. 
Since [<3 n W] = o\uj{o\ + er 2 ) on dcu , we deduce u — 0. From the second line of A*(u,v) — 0, we 
see that Sqv = on 9f2. Since the single layer potential operator SV2 : H _1 / 2 (9f2) h- > H 1 / 2 (i9il) is an 
isomorphism, v = holds. The injectivity of A* is proved. Since 2A = I + C where C is a compact 
operator, we conclude that A has a continuous inverse thanks to the Fredholm alternative. I 

In a similar way, the problem (I57j) is well-posed under some additional assumptions. We define 
the adequate space 



We can state the following result. 

Theorem 8 // we impose the normalizing condition 



v n = I fi 

an Jon 



then there exists one unique couple {{v n )\dui ( v n)\dn) £ H 1 / 2 (9o;) x n]^ 2 (dQ) solution of ([5 
Proof of Theorem [8] Set 



B = 



1 °\ 

-I + fJ,K u ; KqQQu 

2 o-2+o-i 
a x ( 1 

^K duldn ; --I + Kq 

cti + a 2 \ 2 



(61) 



the operator defined on H 1 / 2 (9w) x H^ /2 (dfi). The adjoint B* can be written under the form 



B* 



1 



CTI + (72 



-K, 



o-i ( 1 



Ol + (T 2 



I + K? 



(62) 



In a first step, we begin to show that B* is injective. Let {u,v) € H 1 / 2 (9w) x H 1 ^ 2 ((?r2) be in the 
kernel of B* . We introduce the potential 



Z(x) = ^ — / T(x,y)u(y) + / T(x, y)v(y) ), xG 

o-i + a 2 \J dlAj J dn I 
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We can see that Z is a harmonic function in M. d \(duo U Sf2), satisfying d n Z = on dfl. By the 
properties of the single layer potential, [Z] = Furthermore, a straightforward calculation shows 
that [o~d n Z] = on duo. Hence, Z solves the boundary value problem 

-div (crVZ) = in CI, 
d n Z = on d£l. 

The function is therefore constant in f2. Writing [d n Z] — on duo, we get easily u = and then 
(— 5 + Kq)v — 0. Since the operator A/ — if^ is one to one on H^ 2 (df2), we deduce that v = 0. We 
conclude the proof thanks to the Fredholm alternative. ■ 

Tangential regularity results. Let us consider now the particular case where both a and (3 
are the zero function and (/i,<7i) = (/><?) where / and g are respectively the Dirichlet and Neumann 
boundary data. To recover the tangential regularity of the solution u along duo, we look at the first 
line of (p)6"|) to deduce that 

<7l (Ji 

(«d)iaw = ; Sdnaud n ud\du H ; Kdnduf; (63) 

0"2 + 0"x 01 + 02 

and 

It is easy to deduce that (ud)\du S C 3 ' q (<9lj). Indeed, from (|6"3"|) that we consider as an equation in 
[ud)\du with data / and (d n u d )\ an = g, we see that (/, (d n v d )\dn) belongs to H 1/2 (<9S1) x H~ 1/2 ((9fJ), 
thanks to Theorem [7] 

In order to give a sense to the jump conditions arising in ([6j) , (jTl) , dHJ) , we need to work in space of 
functions of higher regularity. We choose the framework of Holder spaces. We quote [T2] to precise 
the behavior of the layer potentials on these spaces. 

Theorem 9 (Kirsch JW/ ) 

1. If du is of class C 2,a , < a < 1 then the operators and map C^(duj) continuously into 
C 1 ' 13 for all < (3 < a. 

2. Let k £ N with k ^ 0. If duu is of class C k+1,a with < a < 1, then the operators S^, and K M 
map C k 'P(dui) continuously into C k+1,13 (duj) for all < (3 < a. 

3. Let k be an integer. If duo is of class C k+2 ' a , then K* maps C k ^ continuously into C k+1 '^ (du>) 
for all < j3 < a. 

We go back to the proof. Since the two boundaries have no intersection point and since duo is of 
class C 4,Q , it follows that the right hand side of the former equation is of class C 3,a (du>). We then 
conclude the solution of (f6"3"| will be of class C 3 a since the operator 1/21 + /j,K u is an isomorphism 
from C 3 ' a (duj) into itself. With the same arguments, we show straightforwardly that (d n u n )^ duj £ C 2 a . 
About the regularity of the jumps of the second derivative. The equations giving the jump 
conditions [u' d ] and [du' d ] show obviously that [u' d ] and [9 n M^] belong respectively to C 2 ' a (du) and 
C 1,a (dco). Hence, it comes straightforwardly that [u d ] £ C 1,a . With the same arguments, we show 
that [<9„u^] £ C°' Q (see [20] for more details) and then that all the formal computations to get the 
equations describing the second derivative have a sense. 

Remark 1 In a view of a numerical discretization of the state equation, one has to emphasize that 
the choice of a finite elements method seems inappropriate: one should extract tangential derivative 
of high order on the interface duo. The obtained numerical accuracy is not sufficient to incorporate 
the results in an optimization scheme. On the converse, the systems of boundary integral equations 
(|56p and (|57p are well-suited for this kind of computation. Nevertheless, a discussion of adapted 
schemes should be precise and is out of the scope of this manuscript. 
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3.6 Case of Neumann boundary conditions. 

Since the admissible deformation fields have a support with no intersection points with the outer 
boundary, it is a straightforward application of the preceding computations to show that u n solution 
to CHJ)— (J4J) is twice differentiable with respect to the shape. Furthermore, its second order derivative 
u'n belongs to H 1 (il \ W) U H 1 (w) and solves 

Ait" = in lo \ uj and in w, 

[<] = (hi^ n h 2 ^H - hi T .L>nh 2r ) [d n u n ] - (hi, n [d a (u n )' 2 ] + h 2 , n [d n {u n y i \) 
+ (hi r .Vft,2,n + h 2r .V/ii in ) [d n u n ] on dui, 
[ct<9„<] = div T (h 2 ,„ [crV T (u„)i] + hx, n [crV r (u„) 2 ] + hi T .Dn.h 2r )[crV r u„]) 

-div T ((h lr .V T /i 2 ,„ + V r ft,i : „.h 2r + h 2 , n hiji{2Dn - HI)) [<tV t m„]) on du>, 
d n u" — on 5f2; 

(65) 

where we use the notations of Theorem [5J 

4 Second order derivatives for the criterion. 
4.1 Proof of Theorem [3J 

The differentiability of the objective is a direct application of Theorem [5] The computation we make 
here is based on the relation 

D 2 J KV {u)(h ll h 2 )=D(DJ KV {w)h 1 )h 2 -DJ KV (w)Dh 1 h 2 . (66) 

To obtain (fTU|) , we compute in a first step the shape gradient in the direction hi. Then, in a second 
step, we differentiate the obtained expression in the direction of h 2 . In the sequel, we adopt the 
notation v — Ud — u n to obtain concise expressions. 

DJ KV (w)hx = ax J div (jVv| 2 hi) + 2Vf.Vu 1 + a% J div(|Vv| 2 hi) + 2Vv.X7v[ 
= ax(Ax + 2Bx) + a 2 (A 2 + 2B 2 ), 

where 

A x = I div(|Vu| 2 hi) Bx = [ Vv.\7v[ 
A 2 = ( div (|Vi;| 2 hi) B 2 = ( S7v.Vv[. 

J UJ J UJ 

Now, we use the classical formulae to differentiate a domain integral to get 

DAx(uj)h 2 = J div ^div (jVv| 2 hi) h 2 ^ + 2div (\7v.Vv' 2 h x 
= - / div f|Vu+| 2 hi) h 2 , n + 2\7v+.\7(v 2 )+ h X ,n) 
DA 2 (w)h 2 = / div (|Vw"| 2 hi) h 2 , n + 2Vir.V(i> 2 )~ hi, n . 
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The terms DBi, i = 1,2 require more precisions. First, we write 
D5 1 (w)h 2 = / div(v^Vwi)h 2 ) + V^.Vt^ + Vv.V(t>i) 2 , 

7«+.V(«0 + /l 2 ,n + 9n« + (M) + )' 2 + 5 + + d n (v' 2 )+(v[Y 



d n v((Un)[y 2 + ~ [dn{udYi{Un)2 + <9 n (u^ (Un)'i 
dQ. 1 



Note that we used the Green formula twice to keep the symmetry in hi and h2 . We also use the fact 
that the derivatives (ud)'i are harmonic in f2 \ ZU to transform the boundary integral on the exterior 
boundary into an integral on the moving boundary. We obtain 



DBi(w)h a = - f Vu+V(«0 + h 2 ,„ + d n v+ ((Mi) 2 ) + - vd n (((u n )' i y 2 )+ 

J did 

\ (d n (v[) + ((u d y 2 ) + + d n (v' 2 ) + (( Ud )[) + a„K)' 1 ) + (« 2 )+ - d Q ((u n y 2 ) + KV 



1 

Idui 2 

By the same methods, we get 



DB 2 {w)h 2 = f Vv-.V{v[)-h 2 , n + d n v-((v[Y 2 )- + l(d n (v[)-(v' 2 )- + d n (v 2 )-(v[)-) 



We regroup the different terms and after some straightforward computations, we obtain 

2 



D (DJ KV (u)hi) (ui)h 2 = - I div( cr|Vu| 2 hi 

J duj 



+2 f r 



aVv. (hi, n Vv' 2 + h 2tn Vv[ 
o [(u d ) 2 d n v[ + (u d )i<9 n «2 - d n (u n ) 2 v[ - a n (u„)i« 
CT5 n ((w„)i) 2 -cri<9 n w + ((wd)i) 2 ) 



In order to compute D 2 Jjfv(a;)(hi, h 2 ), the first order derivative of the Kohn-Vogelius objective is 
needed. It can be written as follows: 



DJ KV (w)h = - a\Vv\ 

J duj 



h n + 2 



c)uJ 



(Td n (u n )' - aid n v + (u d )' 



Gathering (|66p . (|4Tj) and (|42|) . we write the second derivative of the Kohn-Vogelius criterion as: 



D 2 J KV (uj)(h u h 2 ) 



div 



c)lj 



cj|Vw| 2 hi 



a\Vv\ 



(Dhxh 2 ).n 



a (iu d y 2 d n v'l + (u^dnVl, - d a (u n y 2 v' x - 3 n (Wn)iw2 

aVv. (hi in S7v' 2 + h 2 , n Vv'^j +v ad n (u n )"_ 2 - vid n v A 



Ml 



Let us give a more simplified version for the first term. We decompose the field h 2 into normal vector 
and tangential parts and we use (|43[) . After some elementary computations, we obtain 



div 



cr|Vu| 2 hi 



crlVvl 



(Dhih 2 ).n 



c)uj 



a\Vv\ 2 (hi r .V T /i 2 ,„ +h 2r .V T /ii,„ - h 2T .Unh lT ) - / d n ( a\Vi 

J duj 



hi n h 2n . 
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Finally, the second order derivative of the Kohn-Vogelius objective becomes: 



^ 2 j^M(h l5 h 2 ) 



3oj 



crlVd 



(h lT .V T h 2 , n + h 2T .V T /ii,„ - h 2T .Dnh lT ) 



ctIVwI 



o-V 



CT ( (u<i) 2 9nUi + (Ud)'id n v' 2 - dniUn)^ - 9n(u„)i« 2 



(67) 



2 / u 



0-d a (Un)" 



fad)",: 



4.2 Analysis of stability. Proof of Theorem 

Now, we specify the domain w that is assumed to be a critical shape for Jkv- Moreover, we assume 
that the additional condition Ud — u n holds. To emphasize that we deal with such a special domain, 
we will denote it u>*. The assumptions mean that the measurements are compatible and that lo* is 
a global minimum of the criterion. From the necessary condition of order two at a minimum, the 
shape Hessian is positive at such a point. 

Let us notice that only the normal component of h appears. Let us also emphasize that there 
is no hope to get h = from the structure theorem for second order shape derivative ([BJ). The 
deformation field h appears in D 2 Jxy(^*)(h, h) only thought its normal component h n since u>* is 
a critical point for Jkv- This remark explains why we consider in the statement of Theorem [5] the 
scalar Sobolev space corresponding to the normal components of the deformation field. 

We now prove Theorem [SJ From (|67[) . we deduce 



Djly{w*)[h, h] 



a (u' d dnv' — d n u' n v' 

2 H / I ( u 'd + ~ u « + ) div r [h n V T u d ) -d n ujh n d n (u d - u' n ) 

Jaw \ a i 

2 [a] ( (u'£ - u'+, div T (h n X7 T u d ) ^ - -^-^d n u d h n , d n (u' d - u' n ) 



(68) 



where (,) denotes the duality between H 1 ^ 2 (8lo*) x H 1 l 2 (du!*) . Let us introduce the operators 



Ti : H l / 2 {dw*) 
h 

h 



H-^idw*) Mi : H^ 2 (doj*) -> H 1 /' 2 (dijj*) 
div r {h n V T u d ) h i * u' d + - u'+ 



h n 8r 



tnv n U d 

The Hessian can then be written under the form 



h i-> d a \u 



D 2 J KV (u J *){\\v) =2[a] I/m^T^)) - — /V 2 (h), M 2 (h) 



c 2 \ 

From the classical results of Maz'ya and Shaposhnikova on multipliers ([2], [12]), we get easily that 
T\ and T 2 are continuous operators. In fact, the compactness of the Hessian is a consequence of the 
fact that both operators M\ and M 2 are compact. We use a regularity argument : we remark that 
Mi is the composition of the operators: 



Ri : H 1 / 2 ^*) 
h 



Hl /2 (dn) and R 2 : Hl /2 (dn) 



H 1 / 2 (duj*) 
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where ip is the trace on doj* of \P solution of 



-A* = in fi \ cJ* and in uj* 

[*] = on <9w*, 

[crdn*] = on <9w*, 

* = 6 on c?Sl. 



While R\ is a continuous operator, we prove that R 2 is compact. Let us express u\q u <. 
the integral formula of u to obtain: 



ip. We use 



1 a x 

-I + /i-fC w * ; Sondu* 



C2 + CTl 
<7l 



02 + CTl 



(9 n u)\ 9 Q 



0\ 



(Tl + 02 



Kqq,Ou 

1 



The matricial operator arising in this equation appeared also in 
thanks to Theorem [7] Let us express u\q u * = ip: 



1 _ x 



o-i 



CTl + 0- 2 



56l). It has a continuous inverse 



0. (69) 



Since the operators i^anaw* and SgQg^* are compact, the operator i? 2 is compact, hence Mi is 
compact. The proof of compactness of M 2 is similar. Let us mention that a similar strategy of proof 
can be found in [5]. 

The natural question is then to quantify how is this optimization problem ill-posed. This question 
is directly in related to the rate at which the singular values of the Hessian operator are decreasing. 
Equation (|69|) shows that this rate is the one of the operators Kq^q u * and Sgndu* ■ Now, since for 
every u G H 1 / 2 (917), the functions KgudwU and Sgnawu are harmonic outside of dfl and therefore in 
fi, their restrictions on dui* are as smooth as duj* . We conclude that if dui* is C°° then the restriction 
belongs to each R s (dw*) for s > 1/2 then that if A„ denotes the n th eigenvalue of D 2 Jkv(^*), then 
A„ = o(n~ s ) for all s > 0. 
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